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Cities have been expanding rapidly worldwide, especially over the past few decades. 
Mapping the dynamic expansion of impervious surface in both space and time is 
essential for an improved understanding of the urbanization process, land-cover 
and land-use change, and their impacts on the environment. Landsat and other 
medium-resolution satellites provide the necessary spatial details and temporal 
frequency for mapping impervious surface expansion over the past four decades. 
Since the US Geological Survey opened the historical record of the Landsat 
image archive for free access in 2008, the decades-old bottleneck of data limita- 
tion has gone. Remote-sensing scientists are now rich with data, and the challenge 
is how to make best use of this precious resource. In this article, we develop 
an efficient algorithm to map the continuous expansion of impervious surface 
using a time series of four decades of medium-resolution satellite images. The 
algorithm is based on a supervised classification of the time-series image stack 
using a decision tree. Each imerpervious class represents urbanization starting in 
a different image. The algorithm also allows us to remove inconsistent training 
samples because impervious expansion is not reversible during the study period. 
The objective is to extract a time series of complete and consistent impervious 
surface maps from a corresponding times series of images collected from mul- 
tiple sensors, and with a minimal amount of image preprocessing effort. The 
approach was tested in the lower Yangtze River Delta region, one of the fastest 
urban growth areas in China. Results from nearly four decades of medium-res- 
olution satellite data from the Landsat Multispectral Scanner (MSS), Thematic 
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Mapper (TM), Enhanced Thematic Mapper plus (ETM+) and China-Brazil Earth 
Resources Satellite (CBERS) show a consistent urbanization process that is consis- 
tent with economic development plans and policies. The time-series impervious 
spatial extent maps derived from this study agree well with an existing urban extent 
polygon data set that was previously developed independently. The overall map- 
ping accuracy was estimated at about 92.5% with 3% commission error and 12% 
omission error for the impervious type from all images regardless of image quality 
and initial spatial resolution. 


1. Introduction 

Due to the rapid growth of population and economic development, more people are 
now living in cities than in rural areas for the first time in human history (United 
Nations, 2008). As a result, cities are sprawling rapidly into their surroundings (Seto 
et al. 2002, Schneider et al. 2005, Seto et al. 2010). This is particularly true for China 
because of its large population and the rapid economic growth over the past 30 years 
(Ji et al. 2001, Liu et al. 2005). China’s urban population increased from 20% in 
1982 to 46% in 2010 (Peng 201 1). Urban environments support an improved quality 
of life for most of their inhabitants. One may also see urbanization as advantageous 
to industrial growth and commercial activity. However, excessive expansion of imper- 
vious surface can cause a series of environmental problems, including but not limited 
to loss of habitats for wildlife (McKinney 2002), change of regional climate and water 
cycles (Zhou et al. 2004), loss of fertile land for food production (Seto et al. 2002) 
and altering ecosystem goods and services (Groffman et al. 2003, Grimm et al. 2008). 
These impacts have been recognized and studied for decades (Kalnay and Cai 2003), 
but a complete understanding has been hindered by a lack of consistent geospatial 
data on global urbanization. Clearly, self-consistent and continuous time-series maps 
of impervious surface are needed for this effort over any geographic area concerned. 
A time series of impervious surface maps also provides independent information for 
decision-makers to monitor land-use and land-cover changes, to evaluate the imple- 
mentation of land-use policy (Netzband et al. 2007) and also to potentially model the 
future expansion of cities (Herold et al. 2003, Jantz et al. 2005). 

A characteristic change associated with urbanization is the expansion of impervious 
surface. Satellite remote sensing provides the only viable option to detect and monitor 
impervious surface from space in an efficient, affordable and timely manner. There 
are numerous satellite products that offer the potential to monitor urban growth. 
The night-time lights of the world product derived from the Defense Meteorological 
Satellite Program (DMSP) Operational Linescan System (OLS) was used to produce 
the global urban data set at 1 km x 1 km spatial resolution (Elvidge et al. 1997, 
Imhoff et al. 1997). Moderate Resolution Imaging Spectroradiometer (MODIS) land- 
cover products were combined with DMSP nightlight and population data to produce 
a synergized 1 km x 1 km spatial resolution global urban map (Friedl et al. 2002, 
Schneider et al. 2003). Global urban maps at coarse resolution can cover large areas 
and be updated frequently. However, due to the complex nature of urban landscapes 
and inherent resolution of human activity, global urban maps at the coarse spatial res- 
olution have been difficult to use within many applications and for policy makers at 
local to regional scales (Small 2003). This is also due to their limited mapping accuracy 
and difficulty in accurately describing finer resolution objects such as roads, buildings 
and other infrastructure. 
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Medium-resolution optical satellite imagery possesses unique advantages in map- 
ping urban areas accurately and over time. It provides the detailed information with 
sufficient spatial and temporal resolution for detecting urban expansion. Sensors 
onboard the Landsat series of satellites, including the Multi-Spectral Scanner (MSS), 
Thematic Mapper (TM) and Enhanced Thematic Mapper plus (ETM+), have been 
providing Earth observation data continuously since the early 1970s and form the 
cornerstone for medium-resolution remote sensing (Townshend et al. 1991, Loveland 
and Shaw 1996). Many urbanization mapping activities at local, regional and conti- 
nental scales have been conducted at this spatial resolution (Seto et al. 2000, 2002, 
Jantz et al. 2005, Schneider et al. 2005). As part of the National Land Cover 
Dataset (NLCD), impervious surface maps were produced for the USA from Landsat 
data for the first time (Homer et al. 2004). Data from Landsat have also been 
used successfully in land-use and land-cover change detection (Gutman et al. 2004). 
However, the majority of urban land-cover and land-use studies have tended to 
focus on a single image at one time in the past (Schneider et al. 2003, Wu and 
Murray 2003, Song 2005). Applications using Landsat data are surging thanks to 
the availability of free Landsat data from the US Geological Survey (http://glovis. 
usgs.gov/) since 2008 (Woodcock et al. 2008). Therefore, approaches that take advan- 
tage of the moderate-resolution image time series and produce consistent maps for 
urban expansion will be highly valuable for urban ecologists, city and regional plan- 
ners and remote-sensing scientists monitoring land-cover /land-use changes. These 
data sets would also allow for more detailed studies to determine the impact or 
consequences of such changes on the water, carbon and energy cycles, and even 
climate. 

Continuity of medium-resolution (or Landsat-like) data is critical for monitoring 
land-cover /land-use change. However, the failure of the Scan-Line Corrector (SLC) 
onboard the Landsat 7 satellite in 2003 caused a loss of ~25% of the data toward 
the edges of each image, and Landsat 5 suspended operations in November 2011. 
Although the new Landsat Data Continuity Mission (LDCM) satellite is expected in 
orbit in early 2013, maintaining the continuity of Landsat style medium-resolution 
remotely sensed data is precarious. This situation highlights the need to combine 
the capabilities of existing international sensors to provide a more robust observa- 
tional record. The China-Brazil Earth Resources Satellite-2 (CBERS-2) was launched 
on 21 October 2003. It has a similar satellite orbit (778 km) and equatorial cross- 
ing time (at local time 10:30 am) as Landsat and carries the Charge Coupled 
Device (CCD) and Infrared MultiSpectral Scanner (IRMSS) cameras (Wang et al. 
2000, Sausen 2001). The CCD visible and NIR bandwidths are almost identical to 
Landsat but with a 20 m spatial resolution. It is a valuable source of data to inte- 
grate with the historical Landsat data record and/or to supplement/enhance this 
record. 

In this study, we develop an efficient mapping algorithm to map impervious surface 
for the lower Yangtze River Delta region of China from a time series of images from 
Landsat MSS, TM, ETM+ and more recent CBERS CCD data, spanning almost four- 
decades. Our aim is to produce a time series of consistent impervious surface maps 
from these images for this rapidly urbanizing region of the world, as a precursor to 
a quantitative assessment of regional urbanization. Although the mapping algorithm 
was developed with the lower Yangtze River Delta as a case study, the approach is not 
location specific and can be used elsewhere with a similar set of moderate-resolution 
images. 
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2. Study area and data 

Our research area is located in the lower Yangtze River Delta, China, encompassing 
Suzhou, Wuxi and Changzhou, which are the three major cities within the Landsat 
frame. The Yangtze River Delta is one of the areas with fastest economic growth since 
China opened its door to the world in the late 1970s. Accompanying the rapid eco- 
nomic development and population growth is the dramatic urbanization in the area. 
For example, the non-agrarian population in Suzhou city increased by about 22% from 
1982 to 1985, 25% from 1985 to 1990 and 29% from 1990 to 1995. Industrialization 
had great impacts on small towns in this area too. Resulting from three major eco- 
nomic stimulus policies from the 1980s to 2000s, this area showed clear spatial and 
temporal urbanization patterns associated with each policy change (Ma et al. 2008). 
There were only three cities (Changzhou, Wuxi and Suzhou) in the region in the late 
1970s, each of which had population well under one million. Each of these three cities 
has a population of several million today. In addition, there are five new cities today, 
each with a population above 400 000, and 205 towns, each with a population between 
10 000 and 400 000 (Ma et al 2008). 

To detect the urban expansion in this area, we selected eight Landsat images (WRS- 
2 path = 119 and row = 38; central location: 120.31° E and 31.75° N) from 1973 
(MSS) to 2005 (ETM+) and a CBERS image (provided free from our collaborators) 
in 2006 in this study (table 1). Most of the selected scenes are very clear except for 
the ETM+ image acquired on 26 May 2002. This ETM+ image looks hazy in the 
eastern part of the image. We included this hazy image to examine the robustness of 
the algorithm. We also included a gap-filled SLC-off ETM+ image (2005) in the time- 
series data set to test whether an SLC-off image would cause consistency problems for 
impervious surface mapping. A CBERS-2 CCD image acquired on 3 March 2006 was 
used to examine whether Landsat-like data could be integrated with Landsat data for 
change detection. Due to the proximity in time, the CBERS CCD image could be used 
in place of the Landsat 7 ETM+ SLC-off image to understand whether the SLC-off 
image (acquired a year ago on 3 1 March 2005) could be used for impervious surface 
mapping. 

3. Approach 

The essence of our approach is to take a time series of images from multiple sensors 
over a particular spatial domain as an integrative whole, thus allowing the removal 
of inconsistent changes. Many change detection methods for urban growth analysis 


Table 1. Landsat (Lnd) and CBERS data used in this study (WGS path/row = 119/38). 


Image index 

Date 

Satellite 

Sensor 

Data quality 

1 

16 November 1973 

Lnd 1 

MSS 

Good 

2 

25 May 1979 

Lnd 3 

MSS 

Good 

3 

8 May 1984 

Lnd 4 

MSS 

Good 

4 

23 July 1991 

Lnd 5 

TM 

High 

5 

3 August 1995 

Lnd 5 

TM 

High 

6 

4 May 2000 

Lnd 7 

ETM+ 

High 

7 

26 May 2002 

Lnd 7 

ETM+ 

Some cloud 

8 

31 March 2005 

Lnd 7 

ETM+ 

Gap-filled 

9 

3 March 2006 

CBERS 

CCD 

Good 
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process each image independently and then employ pixel-by -pixel comparison to iden- 
tify land-cover /land-use changes and, thus, may suffer from the problem of error 
propagation in the overlay analysis (Weng and Lu 2009). A schematic outline of our 
processing approach is presented in figure 1 . 

First, remotely sensed data from different sensors were resampled, co-registered and 
orthorectified through the Automated Registration and Orthorectification (AROP) 
package developed by Gao et al. (2009) using an orthorectified Landset image 
acquired on 23 July 1991 as the base image. The AROP package can handle input 
images with different spatial resolutions and/or projections. It combines multiple 
resampling processes (e.g. reprojection, co-registration and orthorectification) into 
one combined process and thus reduces errors from multiple resampling (Gao et al. 
2009). By using a common Landsat image as the base, all input moderate-resolution 
data including Landsat MSS and CBERS CCD data were processed to the same 



Figure 1 . Processing flow chart for mapping continuous impervious surface from time-series 
medium-resolution remote-sensing images. 


Downloaded by [Nasa Goddard Space Flight Center] at 06:47 09 July 2012 


7614 


F. Gao et al. 


grid and spatial resolution. The input data can be in digital number, radiance, 
top-of-atmosphere (TOA) reflectance or atmospherically corrected surface reflectance. 
In this study, we calibrated MSS data to TOA reflectance. The TM and ETM+ images 
were processed to surface reflectance using the Landsat Ecosystem Disturbance 
Adaptive Processing System (LEDAPS) (Vermote et al 2002, Masek et al 2006). 
Because calibration coefficients for CBERS data were not available for this study, we 
used CBERS data in original digital numbers (DNs). All images in the time series were 
stacked into a single multi-temporal composite image before classification. Because 
the training data obtained from an image layer in the stack will only be used for that 
layer during the classification process, the unit of the remotely sensed data would not 
matter, nor would the units of the image layers in the stack need to be consistent as 
long as each layer of the training data and the corresponding layer from the image 
stack to be classified are on a common measurement scale (Song et al 2001). 

In the second step, impervious samples were selected from each image. This was 
done manually for each image based on visual inspection of the brightness of impervi- 
ous surface. Typically, older impervious surfaces are dark in the image, whereas newer 
surfaces are relatively brighter. We assumed that the process of creating impervious 
surface was irreversible during the time of our study (i.e. the impervious surface cannot 
revert to the pervious state). Therefore, the impervious samples from the earlier years 
could be used in the later years. We selected two to eight impervious classes in differ- 
ent brightness levels as training samples from each image. The number of impervious 
classes also depended on the image quality. For example, fewer impervious classes 
could be selected on the MSS images because of their lower signal-to-noise ratio, and 
thus lower brightness contrasts among impervious classes. For each image, we focused 
on new impervious types, as these are not represented in the earlier images. Figure 2 
shows the combined urban training samples on a subset area of the 2006 CBERS CCD 
image covering the city of Wuxi. Polygons delineating urban training samples in the 
earlier image can be used in this image for urban training data collection. 

After collecting a representative set of impervious surface types for all images, we 
computed a normalized spectral distance (similarity) to impervious surface on a pixel- 
by-pixel basis for each image using the impervious surface samples from the same 
image. The smallest distance from all classes is computed using equations (1) and 
(2) and a single minimum distance map was created for each image as: 


where 3c ; and er, are the mean and standard deviations of impervious type i, respec- 
tively; n is the total number of impervious classes within the image; Xj is the pixel value 
for band j in either DN, TOA reflectance or surface reflectance depending on the input 
data; and m is the number of bands. Because the impervious samples included various 
impervious types, such as different urban building materials, surface and illumina- 
tion conditions, the minimum distance map rendered the heterogeneous impervious 
surface types into a simple distance map. 


n 


S m in = min].?,-}, 


/=i 


( 1 ) 



( 2 ) 
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Figure 2. A subset window (the city of Wuxi centred at 31° 33' N, 120° 18' E) showing imper- 
vious samples selected in different years overlaid on the 2006 CBERS CCD image. Note that 
samples for 2002, 2005 and 2006 are too far way to be shown. 



Figure 3. Impervious spectral distance (similarity) map for the Landsat TM image on 3 August 
1995 (the entire study area with same centre as figure 2). The darker pixels in (b) show shorter 
spectral distances to impervious samples and thus are more likely to be an impervious type. The 
three biggest cities in this area from northwest to southeast are Changzhou, Wuxi and Suzhou, 
respectively, (a) Landsat 5 TM, 3 August 1995 and (b) spectral distance to impervious samples. 


Figure 3 shows an example of the spectral distance map for the Landsat 5 TM 
image acquired on 3 August 1995. The darker tone in figure 3(b) shows smaller spec- 
tral distances (i.e. higher similarity) to impervious surface. Only those samples with 
high confidence (either low or high spectral distance) were used as training samples. 
In this study, we considered pixels with Amin < 0.75 (standard deviation) as ‘confident 
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impervious’ samples and with iS Vnin >1.5 (standard deviation) as ‘confident pervious’ 
samples based on our inspection of the distance image. No pixel with a distance value 
between 0.75 and 1.5 was included in the training set. The cloud shadows in the lower 
right corner in figure 3(a) produced small false distance values in figure 3(b). This issue 
was resolved by our subsequent time-series consistency checks. The predefined persis- 
tent pervious samples (from the most recent image) and persistent impervious samples 
(from the earliest image) are used in the sample training regardless of their .S’ m i n values. 

Figure 4 shows three typical samples extracted from the time-series distance images. 
Figure 4(a) represents a persistent pervious pixel selected from the most recent image. 
It shows a large spectral distance at all time steps. The persistent impervious sam- 
ple (figure 4b) selected from the earliest image shows small spectral differences in all 
distance images. The expanded impervious sample (figure 4(c)) computed from each 
individual image shows large spectral distance to impervious samples in the first two 
images, and small distance for the rest of the images in the time series, which implies 
the surface type changed to the impervious type during the interval between the second 
and the third image acquisitions. 

In the third step, the time-series spectral distance maps were then stacked into a 
multiple layer distance image set and were used to expand impervious surface samples 
for classification. Based on the criteria given above, i.e. S m m <0.75 for impervious sur- 
face and .S’ mm >1.5 for the pervious surface, we expanded the training samples from 
the distance image stack. These expanded samples would undergo an additional con- 
sistency check. Based on the irreversibility assumption, persistent impervious samples 
should be selected from the earliest image and should remain impervious in all images, 
and persistent pervious samples should be selected from the most recent image because 
if a pixel was pervious in the most recent image it could not be impervious in any of 
the previous images. For transient impervious surface, once a pixel was converted from 
pervious to impervious, it will remain impervious in the later images. Pixels that vio- 
late these rules would be removed from the training set. Given a set of n images, there 
will be a total of (n + 1) classes, including a persistent impervious class, a persistent 
pervious class and (n - 1) transient impervious classes (table 2). Therefore, we will be 
able to generate a single layer class map from the n- layer distance maps. Due to the 
conservative criteria we applied, we only have class membership for the pixels we are 
highly confident of, thus the expanded training samples can be used for classifying the 
whole image. 


(a) 


(b) 


(c) 



Image index 




Image index 


Image index 


Figure 4. Temporal characteristic of changes in different impervious classes from the time- 
series distance images: (a) persistent pervious surface; (b) persistent impervious surface; and (c) 
expanded impervious surface. 
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Table 2. The complete set of impervious surface classes based on the time of change from 
pervious to impervious in the time-series satellite images (1 = impervious and 0 = pervious. 
T, represents remote-sensing observation at time i). 


Class index 

T i 

T 2 

T 3 


Ti 


T 

1 n 

Description 

0 

1 

1 

1 

1 

1 

1 

1 

Persistent impervious from T \ to T n 

1 

0 

1 

1 

1 

1 

1 

1 

Transient impervious since year T 2 

2 

0 

0 

1 

1 

1 

1 

1 

Transient impervious since year T 3 


0 

0 

0 


1 

1 

1 


i- 1 

0 

0 

0 

0 

1 

1 

1 

Transient impervious since year 7) 


0 

0 

0 


1 

1 

1 


n — 1 

0 

0 

0 

0 

0 

0 

1 

Transient impervious since year T n \ 

n 

0 

0 

0 

0 

0 

0 

0 

Persistent pervious T \-T n 


Once the expanded training samples were generated from the / 2 -layer distance image 
stack, we will no longer need the distance maps to be stacked for classification. In the 
final step, the (n + 1) class map would be used to train a non-parametric decision tree 
classifier (C5.0 by RuleQuest) on the multiple temporal composite image stack with 
all the images arranged in the order of time collected. Unlike traditional maximum 
likelihood classification, which computes class membership probability based on sta- 
tistical parameters from training data assuming normal distribution of the image data, 
the C5.0 decision tree is a rule-based classifier. The classification rules are generated 
from training samples. C5.0 does not make assumptions for data distribution. The 
classifier has been widely used for classifying remotely sensed images, including the 
global land-cover product from MODIS (Friedl et al. 2002). An advantage of a deci- 
sion tree approach is that it can handle missing values, such as gaps in non-overlapping 
regions. C5.0 provides several options, such as adaptive boosting, to improve classifi- 
cation accuracy. In this study, we chose to use the default option. The trained decision 
trees were then applied to the image stack for classification. The classification would 
result in a (/? + 1) class map covering the entire image frame. 

4. Results and analysis 
4.1 Urban expansion in time 

Figure 5 shows a subset of different impervious classes (right column) detected from 
the classification, representing the approximate time period when the pervious-to- 
impervious change occurred. The corresponding satellite imagery (left column) for 
each time period is also provided for comparison. Based on visual inspection, the 
impervious maps on the right column generally capture impervious spatial extent 
well. As our process used the time-series images as a whole for training and classi- 
fication, a single or a minority number of images with less than ideal image quality, 
such as thin cloud/hazy contamination, had a very limited effect on the results. In this 
study, the cloud shadows in the 1995 Landsat 5 TM image (figures 5 (a 5) and (fi5)) did 
not affect our analyses on the impervious maps for the previous (1991) and the next 
(2000) time periods. However, the impervious information in the thick cloudy areas of 
the 1995 imagery could not be detected accurately. Then clouds or hazy areas in the 
2002 Landsat ETM+ (figure 5 (al)) image did not show much impact on the imper- 
vious map. The incomplete image information in the 2002 image was compensated by 
the next (2005) or the previous (2000) clear image, indicating that users might be able 
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to produce a complete impervious map from a multi-temporal composite image con- 
taining multiple partially cloudy images during a relatively short acquisition period 
(we assume no significant change of impervious surface during this period). Gap- 
filled Landsat ETM+ SLC-off scenes may also be used with this approach for the 
same reason. Both starting and ending images are critical in determining the persistent 



G3) m 


Figure 5. Time-series impervious maps (right column) as detected from Landsat or CBERS 
images (left column) centred at 31° 26' N, 120° 28' E. Colour scheme for the right column: 
red = impervious; green = pervious; blue = water. («1 ) Landsat 1 MSS, 16 November 1973; 
(al) Landsat 3 MSS, 25 May 1979; (a3) Landsat 4 MSS, 18 May 1984; (a4) Landsat 5 TM, 
23 July 1991; (a5) Landsat 5 TM, 3 August 1995; ( a6 ) Landsat 7 ETM+, 4 May 2000; 
(a 7) Landsat 7 ETM+, 26 May 2002; (o8) Landsat 7 ETM+, 31 March 2005; and (a9) CBERS 
CCD, 3 March 2006. Time-series impervious maps detected in the years 1973 (M), 1979 (hi), 
1984 (hi), 1991 (M), 1995 (h5), 2000 (h6), 2002 (hi), 2005 (h8) and 2006 (b9). 
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Figure 5. ( Continued .) 


impervious and pervious types. In this case study, we used MSS and CBERS data as 
the starting and ending scenes. Figure 5 shows that results from MSS and CBERS 
images with different spatial resolutions are both reliable and useful and can thus be 
used with confidence in the continuous impervious surface mapping together with 
Landsat TM and ETM+ data. 

The total areas of impervious surface for the common region from all input images 
are shown in figure 6. Results show that the impervious area in the Yangtze River Delta 
area increased slowly before the 1990s, started to accelerate during the 1990s and then 
surged from 2000 to 2006. This can also be observed from the detected impervious 
maps shown in figure 5. Before 1991, cities were just ‘islands’ and expanded slowly 
from town centres Between 1991 and 2000, urban expansion started to accelerate 
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Figure 5. ( Continued 0 



but cities were still well separated. From 2002 onwards, big cities, such as Wuxi and 
Suzhou, were essentially connected by new built-up areas. The pattern of impervi- 
ous spatial expansion seems to follow the diffusion-coalescence model proposed by 
Dietzel et al. (2005), i.e. early-stage urban expansion occurs in isolated patches (dif- 
fusion), with further urban expansion eventually filling the gaps between the urban 
centres, making the previously isolated urban areas connected (coalescence). The tem- 
poral development pattern of impervious surface agrees with the urban development 
plan and the three stages of the national ‘open door’ policy in China. The first stage 
started in the year 1984 when China opened 14 coastal cities to serve as import/export 
ports to foreign countries. The Yangtze River Delta and the Pearl River Delta are the 
two designated special economic zones, where regional developments were limited and 
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Time 

Figure 6. Temporal trajectory of total impervious surface area in the study area. The imper- 
vious surface areas before, during and after the 1990s reflect the overall economic development 
stages in China. 


controlled by the central government. The second stage started in the year 1992 when 
China embraced a socialist market economy. Under this policy, the Yangtze River 
Delta area attracted a large amount of foreign direct investment due to its convenient 
location. Many towns in this area accelerated their development as a result. The third 
stage started in the year 2000. Large chunks of previously unconnected urban areas 
coalesced during 2000 to 2003, as can be clearly seen in figure 5. In this period, the 
industrial economy became a dominant sector in the Chinese economy. The detected 
areal changes of impervious types shown in figure 6 showed similar development 
stages. 

4.2 Urban expansion map evaluation 

Evaluating classification results of historical images is always a challenge as one can- 
not go back in time to collect independent ground truth. Therefore, one’s ability to 
evaluate quality relies on the quality of existing maps that are produced either earlier 
or independently. In this study, we used the existing urban polygon maps produced 
independently by Ma et al. (2008). These maps were produced through visual interpre- 
tation of Landsat images. Although time consuming, visual interpretation remains one 
of the most accurate approaches for image classification as the interpreter can use the 
contextual, topological and morphological information that cannot be explicitly used 
with computer algorithms. Figure 7 shows that the impervious maps produced in this 
study agree well with urban polygons produced by Ma et al. (2008) based on a subset 
of our study area for Suzhou, the largest city in the study area. These urban polygons 
were visually interpreted and digitized from Landsat MSS (4 August 1984), Landsat 
TM (23 July 1991), Landsat ETM (4 May 2000) and IRS-P6 LIS (3 May 2005) images. 
The 1991 and 2000 Landsat images are in fact the same as those used in this study. Ma 
et al. (2008) defined urban polygons as residential and commercial areas, industrial 
and economic development zones and open space (e.g. water and green lands) if it is 
completely surrounded by urban land-use types. Therefore, the definitions of an urban 
polygon and impervious surface do not exactly match. The urban polygons show the 
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Figure 7. Urban polygons (yellow) from a previous study in Suzhou city (centred at 31° 18' N, 
120° 38' E) (Ma et al. 2008) overlaid on the original image (left) compared with the same poly- 
gons overlaid on the impervious maps (right) detected from time-series images in this study, 
(ul) Landsat 4 MSS, 8 May 1984; (a2) Landsat 5 TM, 23 July 1991; (a 3) Landsat 7 ETM+, 
4 May 2000; and Landsat 7 ETM+, 31 March 2005. Polygons overlaid on the impervious maps 
detected on 4 August 1984 (M), 23 July 1991 (hi), 4 May 2000 (b3) and 3 May 2005 (M). 
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extent of the urban/city as a whole object, whereas our impervious maps only show 
impervious/built-up types at the pixel scale, and with potentially more spatial detail. 
Some discrepancies can be found between urban polygons and impervious maps due 
to their different definitions. The 1991 impervious map (figure 6 (b2)) shows smaller 
detected impervious areas than urban polygons. The urban polygons are somewhat 
extended to include new industrial development areas in the suburb of Suzhou City, 
while the 1991 Landsat TM image (figure 6 (a2)) still shows strong vegetation signals 
in these suburb areas. 

Table 3 shows the overall accuracies for both impervious and pervious classes. The 
evaluation samples were selected from each individual image randomly. For each 
image, we selected about 600-700 pixels for accuracy assessment. The urban poly- 
gons in the reference data were used as the basic sampling unit. The impervious types 
include various bright (new) and dark (old) objects or areas. The pervious types are 
mostly located in suburban areas. As samples are purely based on the visual interpre- 
tation of images at the same resolution, the accuracy assessments should not be over 
interpreted. We are more interested in examining the capability of the algorithm and 
tended to select evaluation samples that are hard to distinguish. The overall accuracies 
may be better than what table 3 shows if areas of crops and water are included pro- 
portionally in the accuracy assessment. The overall accuracies vary in different years 
from 87% to 98%. 

A combined classification error matrix is shown in table 4. This includes assessment 
from all input images, including MSS, TM, ETM+ and CBERS. The omission error 
for the impervious type is about 12.2%. The commission error for the impervious type 
is about 2.7%. The larger omission error may partially be due to the sample selection 
approach (polygon-based). A pixel-based sample selection approach would be more 
appropriate for accuracy assessment as evaluation sample polygons usually contain 
some pervious pixels such as rivers and parks. Due to the lack of a time series of high- 
resolution satellite images or aerial photos, we cannot provide as rigorous an accuracy 
assessment in our study as one would prefer. Considering our assessments are based 
on inputs from different spatial resolution images with different image qualities, the 
results are quite promising. We expect the algorithm to perform better if a data set 


Table 3. Overall classification accuracies for different years. 


Year 1973 

1979 1984 

1991 1995 

2000 2002 

2005 2006 

Accuracy (%) 88.0 

87.3 95.4 

94.3 94.7 

98.1 89.9 

90.8 93.7 

Table 4. Combined classification error matrix from all 

images listed in table 1 . 


Classification results 

Producer’s 

accuracy 

Omission 

errors 

Pixels and percent 

Impervious 

Pervious 

Reference data 
Impervious 

5324 

741 

87.8% 

12.2% 

Pervious 

157 

5710 

97.3% 

2.7% 

User’s accuracy 
Commission errors 

97.1% 

2.9% 

88.5% 

11.5% 


92.5% 
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consisted of more consistent images, such as images from a single sensor (e.g. Landsat 
TM) so that there is less noise from variation in spatial resolution and no variation in 
spectral separablility for the same classes across image layers. 

5. Discussion 

In this article, we successfully used remotely sensed data from four sensors with dif- 
ferent spatial resolutions to produce a consistent impervious surface map over several 
decades. We resampled both MSS and CBERS images to 30 m x 30 m spatial resolu- 
tion to match with that of TM/ETM+ data. Results show that Landsat MSS images 
are acceptable in extending the time series of images for mapping impervious surface 
expansion. 

We noted that the algorithm does not address pixels that are partially urbanized. 
Although the 30 m x 30 m spatial resolution impervious product is much more accu- 
rate at local and regional scales compared with the 1 km x 1 km global products 
and is produced at a high temporal resolution from mixed data streams, there are still 
mixed pixels at this spatial resolution. In fact, mixed pixels exist regardless of what 
the spatial resolution is because the natural world features are predominantly contin- 
uous, and typically not organized as regular grids. A possible alternative solution is 
to extract subpixel impervious fractions (Wu and Murray 2003, Song 2005) for the 
analysis. The choice between a discrete classification and a subpixel fraction method 
depends on the scale and application. 

A possible drawback in using 30 m x 30 m spatial resolution data for change detec- 
tion is its sensitivity to the co-registration accuracy for images obtained from different 
sensors or the same sensor at different times. The location shifts due to image-to-image 
registration errors may cause spurious changes. This problem will be particularly seri- 
ous at the edges or with linear impervious features. Image classification based on the 
time-series image stacks requires very accurate location registration at the sub-pixel 
level. Our checking process for consistent impervious areas can preserve high-quality 
time-series training samples and reduce errors from mis-registration or other sources 
of noise. 

Mapping time series of impervious cover at a coarser resolution (e.g. 60 m) may be 
able to produce higher mapping accuracies (less noisy), but a coarser spatial resolution 
will suffer more from the problem of mixed pixels. 

Despite the sensitivity to registration error, the algorithm developed in this study 
that treats multi-temporal images as an integrative whole does offer several signifi- 
cant advantages. First, it enables the capability to remove inconsistent changes. Not 
all impervious pixels remain completely static. Some impervious pixels may appear 
vegetated such as trees planted along newly built roads and parking lots or on roof 
tops. However, these superficial changes do not alter the fundamental property of 
the imperviousness underneath the vegetation, although urban vegetation does pro- 
vide other benefits, such as relief of surface urban heat island effects (Li et al. 2011). 
These types of pixels can be extremely difficult to account for in conventional classifi- 
cations using a single image. The irreversibility assumption for the imperious surface 
suppressed this type of noise. Second, the algorithm is robust to image noise. The 
images we used in this study are from different seasons. The seasonal variation in 
sun angle, phenology and surface wetness would pose a major challenge to a conven- 
tional classification approach. These types of errors have no effect in the algorithm we 
developed here because the training samples were selected based on relative similarity 
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scores, and the classification was performed on a whole set of images as if it were a sin- 
gle multiple band image. Finally, our approach does not require data normalization 
(e.g. conversion to reflectance or other physical units). Data normalization is a huge 
task in conventional land-cover/land-use classification and change detection applica- 
tions (Song et al. 2001). Our algorithm makes normalization unnecessary for the same 
reasons as given above. 

As this is a supervised classification approach, the classification results depend on 
the quality of training samples. Figure 8 shows the combined time-series impervious 
year and the associated confidence map produced from the same C5.0 decision tree 
classifier. The map shows high confidence over pervious classes such as water and crop- 
lands, and lower confidence (<50%, in red) on some impervious areas in the eastern 
suburb of Suzhou. This may be due to the low-quality image (thin cloud/hazy in the 
2002 image) in this area. The classification can be improved by using high-quality 
images or through an iterative selection on impervious and pervious samples on those 
less confident areas. 

Further validation using different data sources such as high-resolution satellite 
images or aerial photos from various years can give a more precise assessment on the 
classification accuracies. It is always a challenge to validate classification results even 
for single image classification. In the future, we plan to focus on other rapidly urban- 
izing areas with more extensive high-resolution satellite images or aerial photos for 
algorithm testing and validation. There are now rich sources of moderate-resolution 
remotely sensed data that allow us to move from more regional to continental, and 
even global scales to map urbanization (Gutman et al. 2008). Significant challenges 
remain for mapping urbanization over large areas, in terms of not only validation but 
also systematically processing data from multiple time periods, sources/instruments 
and seasons. Our algorithm shows a path forward. 




73 79 84 91 95 00 02 05 06 0% 1-50% 51-75% 76-99% 100% 


(Year) 

(a) (b) 

Figure 8. Urban expansion map (a) and the associated confidence maps (b) may be used 
to assess the classification results and to improve time series by selecting more impervious or 
pervious samples in the less confident areas. 






Downloaded by [Nasa Goddard Space Flight Center] at 06:47 09 July 2012 


7626 


F. Gao et al. 


6. Conclusion 

This article developed a new algorithm to map the continuous urban expansion pro- 
cess from time series of satellite images. The algorithm produced urban expansion 
maps of the lower Yangtze River Delta area of China with high accuracy when com- 
pared with an independent data set. The algorithm included remotely sensed data from 
different sources and spatial resolutions, including Landsat MSS, TM/ETM+ and 
data from CBERS. Because the multitemporal images were treated as an integrative 
whole, no noise removal or calibration was needed. The approach effectively integrated 
Landsat TM and ETM+ data with other imagery in monitoring urban expansion and 
showed a strong potential to extend the satellite data record for land-cover /land-use 
change detection. 

Clouds and their shadows were shown to have minor effects on the change anal- 
ysis. The partially cloudy or gap-filled images were also shown to be useful for the 
time-series analysis using our approach, significantly simplifying the image-processing 
efforts. This algorithm can easily be adapted to other remote-sensing data sources 
and extended to other regions. As Landsat data are currently freely available, it is 
now possible to produce regional, continental or global time-series impervious maps 
from historical Landsat or Landsat -like data. Although significant challenges remain 
in effectively and consistently processing a global time series of medium-resolution 
data, time-series impervious cover maps can provide an improved and more detailed 
data source to evaluate the impact of the urbanization process on the environment as 
well as any potential effects on climate change. 


Acknowledgements 

This work was supported by the US Geological Survey (USGS) Landsat Science Team 
project, the NASA Land Cover and Land Use Change (LCLUC) project and the 
Strategic Priority Research Programme of the Chinese Academy of Sciences, Climate 
Change: Carbon Budget and Relevant Issues, Grant No. XDA05050106. 


References 

Dietzel, C., Oguz, H., Hemphill, J.J., Clarke, K.C. and Gzaulis, N., 2005, Diffusion and 
coalescence of Houston Metropolitan Area: evidence supporting a new urban theory. 
Environment and Planning B - Planning and Design , 32, pp. 231-246. 

Elvidge, C.D., Baugh, K.E., Kihn, E.A., Kroehl, H.W. and Davis, E.R., 1997, Mapping 
city lights with nighttime data from the DMSP Operational Linescan System. 
Photogrammetric Engineering and Remote Sensing, 63, pp. 727-734. 

Friedl, M.A., McIver, D.K., Hodges, J.C.F., Zhang, X., Muchoney, D., Strahler, A.H., 
Woodcock, C.E., Gopal, S., Schneider, A., Cooper, A., Baccini, A., Gao, F. 
and Schaaf, C., 2002, Global land cover from MODIS: algorithms and early results. 
Remote Sensing of Environment, 83, pp. 135-148. 

Gao, F., Masek, J. and Wolfe, R., 2009, An automated registration and orthorectification 
package for Landsat and Landsat-like data processing. Journal of Applied Remote 
Sensing, 3, 033515, doi: 10. 1117/1.31 04620. 

Grimm, N.B., Faeth, S.H., Golubiewski, N.E., Redman, C.L., Wu, J., Bai, X. and Briggs, 
J.M., 2008, Global change and the ecology of cities. Science, 319, pp. 756-760. 

Groffman, P.M., Bain, D.J., Band, L.E., Belt, K.T., Brush, G.S., Grove, J.M., Pouyat, 
R.V., Yesilonis, I.C. and Zipperer, W.C., 2003, Down by the riverside: urban riparian 
ecology. Frontiers in Ecology and the Environment , 1, pp. 315-321. 


Downloaded by [Nasa Goddard Space Flight Center] at 06:47 09 July 2012 


Mapping impervious surface with time series satellite images 


7627 


Gutman, G., Byrnes, R., Masek, J., Covington, S., Justice, C., Franks, S. and JJeadley, 
R., 2008, Towards monitoring land-cover and land-use changes at a global scale: the 
global land survey 2005. Photogrammetric Engineering and Remote Sensing , 74, pp. 
6 - 10 . 

Gutman, G, Janetos, A.C., Justice, C.O., Moran, E.F., Mustard, J.F., Rindfuss, R.R., 
Skole, D., Turner II, B.L. and Cochrane, M.A., 2004, Land Change Science: 
Observing, Monitoring and Understanding Trajectories of Change on the Earth ’s Surface, 
457 p. (Dordrecht: Kluwer Academic Publishers). 

Herold, M., Goldstein, N.C. and Clarke, K.C., 2003, The spatialtemporal form of urban 
growth: measurement, analysis and modeling. Remote Sensing of Environment, 86, pp. 
286-302. 

Homer, C., Huang, C., Yang, L., Wylie, B. and Coan, M., 2004, Development of a 
2001 national landcover database for the United States. Photogrammetric Engineering 
and Remote Sensing, 70, pp. 829-840. 

Imhoff, M.L., Lawrence, W.T., Stutzer, D.C. and Elvidge, C.D., 1997, A technique for 
using composite DMSP/OLS ‘city lights’ satellite data to map urban areas. Remote 
Sensing of Environment, 61, pp. 361-370. 

Jantz, P, Goetz, S. and Jantz, C., 2005, Urbanization and the loss of resource lands in the 
Chesapeake Bay Watershed. Environmental Management, 36, pp. 808-825. 

Ji, C.Y., Liu, Q., Sun, D., Wang, S., Lin, P. and Li, X., 2001, Monitoring urban expan- 
sion with remote sensing in China. International Journal of Remote Sensing, 22, pp. 
1441-1455. 

Kalnay, E. and Cai, M., 2003, Impact of urbanization and land-use change on climate. Nature, 
423, pp. 528-531. 

Li, J., Song, C., Cao, L., Zhu, F., Meng, X. and Wu, J., 2011, Impacts of landscape struc- 
ture on surface urban heat islands: a case study of Shanghai, China. Remote Sensing of 
Environment, 115, pp. 3249-3263. 

Liu, J., Tian, H., Liu, M., Zhuang, D., Melillo, J.M. and Zhang, Z., 2005, China’s changing 
landscape during the 1990s: large-scale land transformations estimated with satellite 
data. Geophysical Research Letters, 32, L02405, doi:10.1029/2004GL021649. 

Loveland, T.R. and Shaw, D.M., 1996, Multi-resolution land characterization: building 
collaborative partnerships. In GAP Analysis: A Landscape Approach to Biodiversity 
Planning, J.M. Scott, T.H. Tear and F.W. Davis (Eds.), pp. 79-85 (Bethesda, MD: 
American Society for photogrammetry and Remote Sensing). 

Ma, R., Gu, C., Pu, Y. and Ma, X., 2008, Mining the urban sprawl pattern: a case study on 
Sunan, China. Sensors, 8, 6371-6395. 

Masek, J.G., Vermote, E.F., Saleous, N.E., Wolfe, R., Hall, F.G., Huemmrich, E, Gao, 
F., Kutler, J. and Lim, T.K., 2006, A Landsat surface reflectance data set for North 
America, 1990-2000. IEEE Geoscience and Remote Sensing Letters, 3, pp. 69-72. 

McKinney, M.L., 2002, Urbanization, biodiversity, and conservation. Bioscience, 52, pp. 
883-890. 

Netzband, M., Stefanov, W.L. and Redman, C., 2007, Applied Remote Sensing for Urban 
Planning, Governance and Sustainability (Berlin: Springer). 

Peng, X., 2011, China’s demographic history and future challenges. Science, 333, pp. 
581-587. 

Schneider, A., Friedl, M.A., McIver, D.K. and Woodcock, C.E., 2003, Mapping 
urban areas by fusing multiple sources of coarse resolution remotely sensed data. 
Photogrammetric Engineering and Remote Sensing, 69, pp. 1377-1386. 

Schneider, A., Seto, K.C. and Webster, D.R., 2005, Urban growth in Chengdu, Western 
China: application of remote sensing to assess planning and policy outcomes. 
Environment and Planning & Design, 32, pp. 323-345. 

Seto, K.C., Kaufman, R. and Woodcock, C.E., 2000, Landsat reveals China’s farmland 
reserves, but they’re vanishing fast. Nature, 406, p. 121. 


Downloaded by [Nasa Goddard Space Flight Center] at 06:47 09 July 2012 


7628 


F. Gao et al. 


Seto, K.C., Sanchez-Rodriguez, R. and Fragkias, M., 2010, The new geography of 
contemporary urbanization and the environment. Annual Review of Environment and 
Resources, 35, pp. 167-194. 

Seto, K.C., Woodcock, C.E., Song, C., Huang, X., Kaufmann, R.K. and Lu, J., 2002, 
Measuring land use change with Landsat TM: evidence from Pearl River Delta. 
International Journal of Remote Sensing, 23, pp. 1985-2004. 

Sausen, T.M., 2001, The China-Brazil Earth Resources Satellite (CBERS). International 
Society for Photogrammetry and Remote Sensing (ISPRS) Society, 6, pp. 27-28. 

Small, C., 2003, High spatial resolution spectral mixture analysis of urban reflectance. Remote 
Sensing of Environment, 88, pp. 170-186. 

Song, C., 2005, Spectral mixture analysis for subpixel vegetation fractions: how to incorporate 
endmember variabilities? Remote Sensing of Environment, 95 , pp. 248-263. 

Song, C., Woodcock, C.E., Seto, K.C., Pax-Lenney, M. and Macomber, S.A., 2001, 
Classification and change detection using Landsat TM data: when and how to correct 
atmospheric effects? Remote Sensing of Environment, 75, pp. 230-244. 

Townshend, J., Justice, C., LI, W.., Gurney, C. and McManus, J., 1991, Global land cover 
classification by remote sensing: present capabilities and future possibilities. Remote 
Sensing of Environment, 35, pp. 243-255. 

United Nations, 2008, World Urbanization Prospects: The 2007 Revision (New York: United 
Nations). 

Vermote, E.F., ElSaleous, N. and Justice, C., 2002, Atmospheric correction of the MODIS 
data in the visible to middle infrared: first results. Remote Sensing of Environment, 83, 
pp. 97-111. 

Wang, H., Wang, S., MA, W. and HE, T., 2000, Development and operation of the IRMSS and 
CCD camera for CBERS. In Proceedings SPIE Infrared Technology and Applications 
XXVI, Bjorn F. Andresen, Gabor F. Fulop and Marija Strojnik (Eds.) (Bellingham, 
WA: SPIE), vol. 4130, pp. 9-12. 

Weng, Q. and Lu, D., 2009, Landscape as a continuum: an examination of the urban land- 
scape structures and dynamics of Indianapolis city, 1991-2000. International Journal of 
Remote Sensing, 30 , pp. 2547-2577. 

Woodcock, C.E., Allen, R., Anderson, M., Belward, A., Bindschadler, R., Cohen, 
W., Gao, E, Goward, S.N., Helder, D., Helmer, E., Nemani, R., Oreopoulos, L., 
Schott, I, Thenkabail, PS., Vermote, E.F., Vogelmann, I, Wulder, M.A. and 
Wynne, R., 2008, Free access to Landsat imagery. Science, 320, p. 1011. 

Wu, C. and Murray, A.T., 2003, Estimating impervious surface distribution by spectral 
mixture analysis. Remote Sensing of Environment, 84 , pp. 493-505. 

Zhou, L., Dickinson, R.E., Tian, Y., Fang, I, Li, Q., Kaufmann, R.K., Tucker, C.J. and 
Myneni, R.B., 2004, Evidence for a significant urbanization effect on climate in China, 
Proceedings of the National Academy of Sciences of the United States of America, 101 , 
pp. 9540-9544. 


